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ABSTRACT 

Context. Formation of stellar counter-rotating components in galaxies. 

Aims. In this paper we disentangle two counter-rotating stellar components in NGC 4191 and characterize their physical properties 
such as kinematics, size, morphology, age, metallicity, and abundance ratio. NGC 4191 was selected as candidate to host stellar 
counter-rotating components on the basis of two symmetric peaks in its velocity dispersion two-dimensional map and its irregular 
velocity field. 

Methods. We obtained integral field spectroscopic observations with VIRUS-W and performed a spectroscopic decomposition tech¬ 
nique to separate the contribution of two stellar components to the observed galaxy spectrum across the field of view. We also 
performed a photometric decomposition, modeling the galaxy with a Sersic bulge and two exponential disks of different scale length, 
with the aim of associating these structural components with the kinematic components. We then measured the equivalent width of 
the absorption line indices on the best fit models that represent the kinematic components and compared our measurements to the 
predictions of stellar population models that account also for variable abundance ratio of a elements. 

Results. We have evidence that the line-of-sight velocity distributions (LOSVDs) are bimodal and asymmetric, consistent with the 
presence of two distinct kinematic components. The combined information of the intensity of the peaks of the LOSVDs and the pho¬ 
tometric decomposition allows us to associate the Sersic bulge and the outer disk with the main kinematic component, and to associate 
the inner disk with the secondary kinematic component. We find that the two kinematic stellar components counter-rotate with respect 
to each other. The main component is the most luminous and massive, and it rotates slower than the secondary component, which 
rotates along the same direction as the ionized gas. The study of the stellar populations reveals that the two kinematic components 
have the same solar metallicity and sub-solar abundance ratio, without the presence of significant radial gradients. On the other hand, 
their ages show strong negative gradients and the possible indication that the secondary component is the youngest. We interpret our 
results in light of recent cosmological simulations and suggest gas accretion along two filaments as the formation mechanism of the 
stellar counter-rotating components in NGC 4191. 

Key words, galaxies: abundances - galaxies: kinematics and dynamics - galaxies: formation - galaxies: stellar content - galaxies: 
individual: NGC 4191 


1. Introduction 

The phenomenon of counter-rotation, i.e. the presence of mul¬ 
tiple kinematic components that are rotating in opposite direc¬ 
tions, has been detect ed in a numbe r of galaxies of all mor¬ 
phological types (see ICorsinil l2014l for a review). Counter¬ 
rotating galaxies are classed according to the nature of their 
counter-rotati ng components, i.e., gas vs. star s (e.g., NGC 4546 , 
,Ganettajl987E stars vs. stars (e.g., NGC 4550. lRubin et ajj jj92 : 


Rix et al 


1992ll . and gas vs. gas (e.g., NGC 7332, Fisher et ^ 


I 994 I 1 . For the cases of stellar counter-rotation, a further clas- 


counter-ro tation is visible only in the innermost regions (e.g., 
NGC 448. lKrainovic et al.ll201 ill . 

Thanks to the adv ent of integral field spectroscopic 
surveys like SA.URON (Ide Zeeuw_et_alj _ 2()()2|l^. ATLAS 3D 


sification can be done by looking at the sizes of the decoupled 
structures. There are galaxies with two counter-rotating stellar 
disks of similar sizes (e.g., NGC 4550) and galaxies where the 
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( Cappeljari^t ah 201 ih . CALIFA (Sanchez et al. 201^ . SAMI 
( Bryant et al.ll2015l) . and MANGA (iBundv et al.ll20 15*1) 7 the cen¬ 
sus of counter-rotating galaxies has increased. Indeed, these 
surveys allowed to identify candidate galaxies to host counter¬ 
rotating stellar disks by looking for the kinematic signature given 
by two off-center and symmetric peaks in the stellar velocity 
dispersion in combination with zero velocity rotation measured 
along the galaxy major axis. These kinematic features are ob¬ 
served in the radial range where the two counter-rotating com¬ 
ponents hav e roughly the same luminosity and the i r LOSVDs are 
unresolved (iRix et al. j_992 Bertola et al.l 119921 : IVergani et al.l 
l2007h . Recently. iRrainovic et al.l ( 2011 ) have found 11 galaxies 
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(including NGC 4191) with a double-peaked velocity dispersion 
in the volume-limited sample of 260 nearby early-type galaxies 
gathered by the ATLAS-3D project. 

The current paradigm that explains stellar counter-rotation 
is a retrograde acquisition of external gas and subsequent 
star formatio n dThakar & RvdenI 119961 : iPizzella et al.l 120041 : 
lAlgorrv et al.l 12014 ). Alternative scenarios such as the assem- 
bly of the counter-rotating stellar component from mergers 

or internal in- 


Crocker et al 


K onent 

200 ^ 


jPuerari & Pfenniger 2001 
stabilities (lEvans & Collett 11994 ) have been proposed, too. The 
relatively small number of studied cases favors external origin, 
but it does not yet allow us to distinguish between ga s accretion 
or merging (ICoccato et al.ll2013t IPizzella et al.ll20 A . 

In order to characterize the physical properties of the two 
stellar counter-rotating components and therefore to constrain 
their formation mechanism, we need to disentangle their con¬ 
tribution to the obseryed galaxy spectrum. To this aim, we in¬ 
troduced a spectral decomposition technique that allows one to 
measure the kinematics and p roperties of the stella r populations 
of the decoupled com ponents (ICocc ato et alJ l201 Ih . 

Other parametric (Johnston et al.ll2013l) and non parametric 
dKatkoy & ChilingariMl^l ih techniques haye also been pro¬ 
posed and deyeloped by other teams. All these studies show that 
the stellar component that is rotating along the same direction 
as the ionized gas is younger, less massiye, and has different 
metallicity and abundance ratio with respect to the main galaxy 
component. These results are consistent with the gas accretion 
scenario followed by star formation. We also found that, at least 
for one case, the small counter-rotation obseryed in the inner re¬ 
gions is the “tip of the iceberg” of a much larger structure, which 
is less luminous than the main stellar component and whose real 
exten t can be reyealed on ly by a spectral decomposition (NGC 
3593. IC^to et al.ll2013h . 

In this work, we inyestigate the SO galaxy NGC 4191, 
which is an isolated s ystem located at a distance of 42.4 Mpc 
dTheureau et al.l l2007h . In contrast to the galaxies we studied 
in our preyious works dCoccato et al.ll20lTll2013l: IPizzella et alJ 
l2014h . NGC 4191 was not already known to host two large-scale 
counter-rotating stellar disks. We selected NGC 4191 because 
its stellar yelocity field reyeals an irregular structure with little 
rotation, and its stellar yelocity dispersion field shows two sym¬ 
metric peaks located on opposite sides at ~ 10” from the galax y 
center along the photometric major axis dKrainoyic et al.ll201 111 . 
The latter feature is one of the key obseryables proposed for de¬ 
tecting candidate galaxies with stellar counter-rotation. Contrar- 
ily to its complex kinematic structure, the photometric profile of 
NGC 4191 does not reyeal prominent sub-structures. The light 
distribution of NGC 41 91 has been modeled w ith a single Sersic 
law with index n - 2.4 dKrainoyic et al.ll20r^ . 

The aim of the paper is to study the structure of NGC 4191 
and to inyestigate the presence of multiple kinematic stellar com¬ 
ponents. In particular, we want to test whether NGC 4191 does 
host counter-rotating stellar components and, if present, we want 
to characterize their kinematics and stellar population proper¬ 
ties. Also, we want to inyestigate whether the kinematic com¬ 
ponents can be associated with any photometric sub-structure, 
which might haye been oyerlooked in preyious studies. More- 
oyer, if two counter-rotating components are indeed detected, 
we will haye strengthened the case that the presence of a double- 
peaked signature in the stellar yelocity dispersion field combined 
with the absence of regular rotation can be used as selection cri¬ 
teria for counter-rotating galaxies. 

This paper is structured as follows; Section describes 
the spectroscopic obseryation and the data reduction; Section [3 



Fig. 1. Location of the observed stars (black dots) in the [MgFe]’, Hy0 
(left panel) and the Mg b, (Fe) (right panel) parameter space. The line 
grids give the predi ctions from simple stellar population models by 
iThomas et al.l l l2003h . 


presents our analysis and our results. Finally, we summarize our 
conclusions in Section |4] 


2. Observations and data reduction. 

2.1. Observations 

The observations of NGC 4191 were carried on 1-3 and 6-7 April 
2013 using the VIRUS-W Integral Field Unit (IFU) Spectro¬ 
graph dFabricius et al.ll2()T^ at the 2.7-m Harlan J. Smith Tele¬ 
scope of the McDonald Observatory (Texas, US). The IFU con¬ 
sists of 267 fibres arranged in a rectangular grid, which covers a 
field of view 105" x 55" with a fill factor of 1/3. Each fibre has 
a diameter of 372 on the sky. We used the lower resolution mode 
of the instrument which covers the 4340A - 6040A wavelength 
range with a spectral resolution of 1.57 A FWHM (ixinstr ~ 38 

kms“')[3. 

The observations were dithered to fill the entire field of view. 
We took 6 exposures of 1800 sec in each of the three dither posi¬ 
tions and bracketed and interleaved them with 300 sec sky nods. 
Each of the science exposures was split in two for cosmic ray 
rejection. The total on-target exposure time per fibre is 3 h. We 
recorded bias frames and Hgn-Ne arc lamp exposures for wave¬ 
length calibration on the evenings and mornings before and after 
the observations. We also recorded dome flat exposures to trace 
the fibre positions on the detector and to compensate for fibre-to- 
fibre throughput variation. A set of spectroscopic standard stars 
were also observed with the same instrumental set-up. 

2.2. Data reduction 

T he data reduction follo ws the same methodology as described 
in iFabricius et al.l(l2014l) . In short, after subtraction of the master 
bias the fiber positions are traced by searching along the peaks 
of the spectra in the dome flat frames. Then the arc peak posi¬ 
tions are identified in the master arc frames by searching along 
those traces. The fiber/wavelength to x/y pixel position mapping 
is then modeled and fit — using a standard least squares mini¬ 
mization routine — with a two-dimensional 7-th degree Cheby- 
shev polynomial. We use 28 Hg and Ne spectral lines for the 
fit. 

After the wavelength calibration we average the two cosmic 
ray split science exposures as well as the bracketing sky frames 

’ The spectral resolution of VIRUS-W below 4800 A deteriorates up 
to (Tinstr ^ 70 km s“'. We tested that restricting the wavelength range to 
4800A - 6040A has negligible effects on our results. 
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while rejecting spuriou s events using the PyCosmic routine 
(iHusemann et al.ll2012h . We extract science, sky and dome-flat 
spectra by stepping along the previously determined trace posi¬ 
tions. For each step, the counts in the CCD pixels are assigned 
to spectral pixels according to the areal overlap of the CCD 
pixel with the 7 pixel wide and 0.52 A long extraction aper¬ 
ture. We correct for flber-to-flber throughput variations using the 
extracted dome-flat spectra and remove the sky background by 
scaling the sky spectra by the relative exposure time. 

Then, we combine all flat-fielded and background-removed 
science spectra into one datacube. We impose a regular grid with 
1'.'6 X 1'.'6 large pixels and calculate the areal overlap between all 
circular fiber apertures and the pixels. The flux in each fiber is 
then distributed to pixels according to the overlap. 

Finally, we spatially rebin the final datacube to increase the 
observed signal-to-noise r atio (S /N) using the imple mentation 
for Voronoi Tesselation bv ICappellari & CopinI (l2003h . 


3. Analysis: kinematics, photometry, and stellar 
populations 

In this Section, we measure the stellar and ionized gas kine¬ 
matics of NGC 4191 and its structural components. The stel¬ 
lar library we use is presented in Section 13.11 We first analyze 
the spectra with a non-parametric fit to recover the shape of 
the line-of-sight velocity distribution (LOSVD) in each spatial 
bin (Section I3.21 i. Motivated by the bimodality of the recov¬ 
ered LOSVDs, we apply a spectral decomposition technique to 
investigate and characterize the presence of two kinematically 
distinct stellar components (Section [3.31) . In Sect. 13.3[ we also 
study the structural components of NGC 4191 and investigate 
the connection between the photometric and kinematics compo¬ 
nents. We then measure the properties of stellar populations of 
the two kinematic components; line strength indices, age, metal- 
licity, and abundance ratios of a elements (Section iTAl) . The dis¬ 
tribution of the ionized gas is discussed in Section [33] 




V [km s '] 


Q 

> 

if) 

O 


0.12 

0.10 

0.08 

0.06 

0.04 

0.02 

0.00 



2000 2500 3000 3500 

V [km s”^] 


-40 -20 0 20 40 

V stors [km s"'] 



0.12 


.0.25 

0.10 



0.20 

Q 0.08 



■ ■'> 0.15 

(/I 0.06 

rF 

- 

■ if) 

y 0.04 

1 


g 0.10 

0.02 

J 


0.05 

0.00 


4,^_ 

0.00 



2000 2500 3000 3500 2000 2500 3000 3500 

V [km s”^] V [km s”^] 


Fig. 2. Best fit kinematics for the single stellar component model. 
The two-dimensional maps show the velocity and velocity dispersion 
fields. The plots on the side show some example of reconstructed non- 
parametric LOSVDs in the location of the field of view indicated by the 
dotted lines. The LOSVDs of central regions clearly show a bimodal 
distribution; the LOSVDs in the outer regions show high asymmetry, 
which can be due to the presence of an unresolved secondary kinematic 
component. 


3.1. The stellar library 

The stellar library used in the spectroscopic fit consists of 12 
stars observed with Virus-W using the same observational set-up 
of the galaxy observations. This ensures that the stellar templates 
had the same line spread function of the galaxy, minimizing sys¬ 
tematic effects due to mismatch between the resolution of the 
stars and the galaxy spectra. This is particularly important for 
systems with low velocity dispersion, where the asymmetries in 
the LOSVDs are comparable with asymmetries in the line spread 
function. 

The selected stars covered a large range of the Hy6, Mg b, 
[MgFe]’, and (Fe) parameters space, as shown in Figure [T] 

3.2. Shape of the line-of-sight velocity distributions 

We fitted the spectra with a non-parametric approach in order 
to recover the shape of the LOSVDs and investigate the pres¬ 
ence of multip le kinematic compone nts. The fitting procedure 
is described in iFabricius et akl ( |2014|) , and it extends t he ma xi- 
mum penalized likelihood method by iGebhardt et al.l (l2000l) to 
treat nebular emissions. After fitting and removing the contin¬ 
uum with a 3-rd degree polynomial, we fitted the wavelength 
range from 4905 A- 5325 A. We modeled the LOSVD with 31 
independent velocity channels and linearly combined the set of 
12 stellar spectra from a library (Sect. ITTI l. We successively re¬ 


moved stellar spectra that received small weights. The final fit 
was carried out with six templates. Due to the high S/N of the 
data after binning (typically 90 per A) we did not employ any 
penalization. For a number of bins the derived LOSVDs show a 
clear bimodal structure expressed as a double peak. 

These bifurcated LOSVDs are more evident in the inner re¬ 
gions of the galaxies (r < 10") and are spatially consistent with 
the regions where the sigma-peaks were detected. For the re¬ 
gions outside 10" the LOSVDs are strongly asymmetric, consis¬ 
tent with the interpretation that the galaxy contains two counter¬ 
rotating stellar components. 

Figure |2] shows the shape of the LOSVDs in some spatial 
bins. In order to measure the mean kinematics along the line- 
of-sight, we re-fitted the spectra parametrizing the LOSVD in 
each bin with Gaussian function plus high-order Gauss-Hermite 
moments. The fit was done using a modified version of the Pe¬ 
nalized Pixel Fitting code of lCaPDellari & EmsellemI (l2004h that 
included simultaneous fits to the ionized gas emission lines. Fig. 
|2]shows the velocity and velocity dispersion fields obtained with 
the parametric fit. The two peaks in the velocity dispersion pro¬ 
file discoveredbv iKrainovic et ^(1201 ll) are visible, slightly off- 
setted by ai 1" towards South East with respect to the galaxy 
photometric major axis. The velocity field is complex and it is 
not consistent with a single rotating component, even if high or¬ 
der moments are taken into account in the LOSVD parametriza- 
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tion, providing further support to the presence of kinematically 
decoupled components. 


3.3. Kinematics and photometry of the counter-rotating 
steiiar components 

Driven by the results of Section !?^ we apply a spectral decom¬ 
position technique to investigate and characterize the presence 
of two kinematically distinct stellar components. 

The spectral decomposition, described in Section 13.3.11 re¬ 
turns for each spatial bin the kinematics, the spectra, and the 
normalized flux contribution of the two stellar components. We 
also fit the intensities of the emission lines, and their mean ve¬ 
locity and velocity dispersion. We used the same spectral library 
for both the components, and let the code to select the most ap¬ 
propriate stellar template for each of them. 

In Section 13.3.21 we associate these two kinematic compo¬ 
nents with the galaxy structural components derived from a pho¬ 
tometric decomposition. 


3.3.1. Spectral decomposition technique 


The parametric fit of the stellar kinematics and the separation of 
two kinematic stellar components is done using the spectr al de¬ 
composition technique developed in ICoccato et al.l (|201 Ih . im- 
plemented as modification to the Penalized Pixel Fitting code 
(ICappellari & Emselleml[2004ll . 

The technique builds two optimal templates by linear com¬ 
bination of stellar spectra from an input library (see Sec¬ 
tion I3.1I) . The two templates are convolved with two indepen¬ 
dent LOSVDs, which are parametrized by two independent 
Gaussian functions. The convolved templates are then multi¬ 
plied by Legendre polynomials to account for the shape of the 
continuum. Our implementation includes Gaussian functions 
that are added to fit the ionized gas emission lines Hj3 and 
[O III] AA4959, 5007A. The best fit model is recovered using the 
Levenberg-Mardquart minimization. The procedure allows 
the flux of the two components either to be free parameters in 
the fit or to be constrained to specified values or within a certain 
range. The code returns the best fit values for the stellar kinemat¬ 
ics and the best fit spectra of the two stellar components, which 
will be used in Section I3.4I to derive their stellar populations. 
Also, the code returns the light contribution of one component 
relative to the total galaxy spectrum. If we label the two com¬ 
ponents with “A” and “B”, the code returns Fa - where 

Ia and Ib are the flux of the two components. By construction, 
Fb^I- Fa. 


Our technique has been successfully applie d to disentan¬ 
gle counter-rotating stellar dis ks in galaxies (ICoccato et al.l 
|201 iL |20T^ Pizzella et al.ll2014l) . stars from the bulge and disk 
jFabricius et al.ll2014 ). and stars from the host galaxy and the 
surrounding polar disk (ICoccato et al.ll2014l) . 


For the special case of NGC 4191, in order to obtain a rea¬ 
sonable and stable fit, it is necessary to constrain the flux ratio 
of the two kinematic components by performing an independent 
photometric decomposition (Section r3.3.2l i. This is due to a de¬ 
generacy between the parameters that describe the kinematics 
and the stellar populations o f the two components. The sam e ap- 
proach has been ado pted in iKatkov & ChilingarianI (1201 111 and 
ICoccato et al.l (l2014t) . 


3.3.2. Photometric contraints 

The aim of this section is to perform a photometric decomposi¬ 
tion of NGC 4191 in order to constrain the flux parameter Fa in 
the spectral decomposition code. 

We acquired a g-band image (the closest match to the 
VIRUS-W wavelength range) from the SDSS DR 7 archive an d 
analyzed it using the 2D image-fitting code IMFIT (lErwi nl2015h . 
The fitting procedure includes the convolution for a Moffat point 
spread function (PSF) computed using the median FWHM and 
P values measured from stars in the same image. 

The best fit model is obtained adopting three components; 
a bulge with Sersic profile, one inner exponential disk, and one 
outer exponential disk. In the fitting procedure, the center has 
been kept fixed for all the 3 components. The position angle and 
the ellipticity of the three components are independent, but con¬ 
stant with radius. 

The structural parameters of the best fitting components are 
given in Table[Tl and the best fit model is shown in Figures|3]and 
|4] The residual map of Fig.|3]highlights a spiral-arm pattern, in¬ 
dicating the presence of faint spiral arms and dust in this galaxy, 
despite of its classification as SO. 

Our photometric decomposition represents an improvement 
of th e previous single component model by iKrainovic et al.l 
1201 3l One way to evaluate the relative goodness-of-fit of dif¬ 
ferent models to the same data is via information-theoretic 
measures such as the Akaike Information Criterion (AIC) or 
Bayesian Information Criterion (BIC), both of which are output 
by Imfit at the end of the fitting process. Unlike the statis¬ 
tic, AIC and BIC do not require that models being compared be 
nested. The absolute values of these statistics are irrelevant; what 
matters is the relative value, with lower values indicating a better 
match to the data. Traditionally, A AIC = AICi - AIC 2 < -10 
is taken as an indication for the clear superiority of model 1 
over model 2. We find that the Sersic + two exponential model 
is clearly favored. In the case of the g-band image, the three- 
component model has AAIC = -1442 compared to the Sersic H- 
single-exponential model and AAIC = -15006 compared to the 
single-Sersic model. The values for the BIC are similarly large: 
ABIC = -1398 and -14919, respectively. 

Figure |5] provides additional support to the choice of three- 
components model over one- or two- components models. It 
shows the radial profiles of ellipticity and position angle as mea¬ 
sured from the g-band image and from the various models. We 
can see that the single-Sersic model does a particularly poor job 
of reproducing the position angle radial variation, and the two- 
components model does not reproduce the change of ellipticity 
as nicely as the three-components model for 25" < a < 30" 
and a > 45". For use in constraining the flux parameter F in the 
spectral decomposition, we therefore use the three-components 
model. 

We define “main kinematic component” the sum of the bulge 
and outer disk and “secondary kinematic component” the inner 
disk. The main kinematic components contains the ~ 64% of 
the total galaxy luminosity, whereas the secondary kinematic 
component contains the ~ 36% of the total galaxy luminosity. 
The particular association of the inner disk component with the 
counter rotating kinematic component is driven by its relative 
local luminosity. From the observed peaks in the LOSVDs we 
know that the two components have flux ratios between 1:3 and 
1:1. Other combinations, such as defining the secondary com¬ 
ponent as the sum of the two disks, would have set one of the 
components below a detection limit {F < 0.15, see Appendix 

which is inconsistent with the observations. 
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Table 1. Best fit parameters of the photrometric decomposition of NGC 
4191. 




Sersic bulge 


PA 

13.3 

± 

0.7 

[degrees] 

ell 

0.128 

± 

0.003 


n 

2.87 

± 

0.06 


Pe 

19.7 

± 

0.05 

[mag arcsec^^] 

fe 

111 

± 

0.07 

[arcsec] 

mxoT 

14.89 

± 

0.05 

[mag] 


Inner exponential disk 

PA 

7.8 

± 

0.1 

[degrees] 

ell 

0.384 

± 

0.003 


AO 

19.78 

± 

0.02 

[mag arcsec^^] 

h 

6.54 

± 

0.05 

[arcsec] 

mroT 

14.23 

± 

0.02 

[mag] 


Outer exponential disk 

PA 

-2.7 

± 

0.7 

[degrees] 

ell 

0.189 

± 

0.004 


AO 

2216 

± 

0.01 

[mag arcsec^^] 

h 

20.9 

± 

0.1 

[arcsec] 

mxoT 

14.02 

± 

0.01 

[mag] 



40 20 0 -20 -40 40 20 0 -20 -40 40 20 0 -20 -40 

ARA [arcsec] ARA [arcsec] ARA [orcsec] 


0 40 80 0 40 80 - 0 . 6 - 0.3 0.0 0.3 0.6 

[ADU] [ADU] [ADU] 

Fig. 3. Result of the photometric decomposition of NGC 4191. Left 
panel: SDSS g image; central panel: best fit model; right panel: normal¬ 
ized residuals (data - model)/model. North is up and East is left. 


At each position {x, y) in the field of view, we can define the 
fraction of galaxy surface brightness contributed by the main and 
secondary components as: 



Fig. 4. Results of the photometric decomposition. Filled circles: surface 
brightness radial profile measured on the SDSS g-band image with the 
iraf task ellipse. The blue curve is the best fitting inner disk component, 
the red-dashed curve is the best fitting Sersic bulge, and the red dot- 
dashed curve is the outer disk component. The green curve is the sum 
of all the best fitting components. Red and blue curves represent the 
profiles of the adopted definition for the main and secondary kinematics 
components, respectively. 



0 20 40 60 

Semi major axis [arcsec] 


Fig. 5. Comparison between the radial profile of ellipticity (top panel) 
and position angle (bottom panel) obtained with different photomet¬ 
ric models. Filled circles: surface brightness radial profile measured on 
the SDSS g-band image with the iraf task ellipse. Green line: three- 
components model; orange dashed line: two-components model; ma¬ 
genta dashed line: single-component model. 


Amain (.^^y) 


Asecond(.^? j) 


_ fmain(.^? T) 

fmain(.^? T) "b Aecond(.^. T) 

fbulge(it'r y) "b ^outer disk(x,y) 

^bulge(-^; j) + ^outer disk(jt,y) + Anner disk(-^,j) 
^second(-^9 y) 

Fb = - 

y) ^second(-^; j) 

dinner disk(x,y) 

^outer disk(.*,y) + Anner disk(.*,y) 

1 Amain(i’tr y) 


( 1 ) 

( 2 ) 


Amain + 0.05; Amain IS giveu by: 

p _ 2 ^main(i’t; y) 

H Immnix, y) + 2 /second(jt, y) 

where the sum is done over all the spaxels contributing to a 
given spatial bin. 

3.3.3. Kinematics of the counter-rotating steiiar 
components 


In Figure|6]we show the radial trend along the major axis of 

Asecond ■ 

The results of the photometric decomposition can be now 
used to constrain the spectral decomposition analysis. For a 
given spatial bin, we constrain the parameter Fa in the spectral 
decomposition code to be within the range Amain - 0.05 < Aa < 


Figure |7] shows the kinematics of the two stellar components, as 
obtained from the spectral decomposition fit. 

The main galaxy component, which consists in the combi¬ 
nation of the Sersic bulge and the large-scale exponential disk, 
is more luminous. It is characterized by a rotation amplitude of 
~ 70 km s“*, with some scatter, and an average velocity disper¬ 
sion of ~ 120 km s“*. 
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Semi mojor oxis [arcsec] 

Fig. 6. Radial profiles of the relative flux contribution of the secondary 
stellar components Fsecondary The contribution of the main component 
is Fmain = 1 “ ^secondary by Construction. 


The secondary galaxy component, which is represented by 
the inner exponential disk, is less luminous and it counter-rotates 
with respect the main stellar component. It has a rotation ampli¬ 
tude of ~ 130 km s“' and an average velocity dispersion of ~ 90 
km s“'. The direction of rotation is aligned with that of the ion¬ 
ized gas. 

As additional test, we have performed the spectral decom¬ 
position using the constrains from the two-components photo¬ 
metric model. Although the main kinematic results are consis¬ 
tent with those presented above, the two-dimensional kinematic 
maps obtained with the two-components photometric model 
have higher noise. This is due to the poorer quality of the photo¬ 
metric fit, which translates into larger uncertainty in constraining 
the light content of the counter-rotating stellar components. 

3.4. Stellar populations of the counter-rotating stellar 
components 

We measured the H/l, Mg b, Fe5270, and Fe5335 line-strength in¬ 
dices on the optimal template of each stellar component returned 
by the spectral decomposition code. The spectra are broadened 
with a Gaussian function to match the spectral resolution of the 
Lick system (8.4 A. I Worthev & OttavianilflOQ^ . Therefore, we 
obtained two sets of indices in each spatial bin, one set for each 
kinematic component. 

Measurements were scaled to the Lick system by comparing 
the line strength indices measured in the spectra Lick standard 
stars HD102870 a nd HD 125560 obtain ed with VIRUS-W with 
those published bv l Worthev et alJ(ll994h . The offsets are -0.35 + 
0.13 A (H;6), -0.05 + 0.01 A (Mg b), -0.02 + 0.15 A (Fe5270), 
and -0.06 + 0.12 A (Fe5270). We applied this correction to the 

and Mg b indices, but not for the metal indices, because they 
are negligible considering the measurements errors. 

In Figure 0 we show the combined indices H/3, Mg b 
[MgFe]’, and (Fe) of the two stellar components. The scatter 
in the plots is large, reflecting the uncertainties in our measure¬ 
ments. In particular, we note that the scatter in the measurements 
in the Mg b (Fe) plane of Fig.[8]is larger for the secondary com¬ 
ponent. This is probably due to a combined effect of i) radial 
gradients that are more pronounced for the secondary compo¬ 
nent (see below); and ii) secondary component being on average 
the less luminous, and therefore the errors in the measurements 
are larger. 

Fig. [8] also compares our measurements to th e predictions 
of stellar population models (iThomas et al.ll20()^ . This allows 
us to derive their two-dimensional maps of age, metallicity, and 


abundance ratios of a elements, which are shown in Figure |3 
Figs.[3and|9]do not show the results on those spatial bins where 
the S /N is lower than 25 per A. In Fig. 0we also compute median 
values within elliptical anuli to highlight the presence of radial 
gradients in the properties of the stellar populations. The error in 
each elliptical bin is computed as the maximum between i) the 
scatter of the measurements within that elliptical bin, and ii) the 

weighted mean error 1 / -^2/ i/crj, where cr,- are the individual 
measurements errors in age, [Z/H], and [a/Fe]. 

The properties of the stellar populations of the two counter¬ 
rotating stellar components are very similar, as shown by Figs.|8] 
and|3 Because the scatter in the two-dimensional maps is large 
and a direct comparison is difficult, we compare average quan¬ 
tities. Metallicities and abundances of a elements are consistent 
with solar, without significant radial gradients. Their mean val¬ 
ues are: <[Z/H])n,am = 0.02 + 0.09 dex, <[Z/H])second = -0.1 +0.2 
dex, <[a'/Fe])n,ain = -0.05 + 0.05 dex, and <[a/Fe])second = 
-0.02 + 0.05 dex. On the other hand, the age radial profiles 
show a mild indication of radial gradients. The age of the stars 
in the main component ranges from ~ 12 Gyr in the central re¬ 
gions down to ~ 6 Gyr at 21". Outside 21" the age rises again, 
but the uncertainties at these radii become too large to derive 
conclusive results. One possible explanation is that, if the ob¬ 
served raise in age is real, it might be due to accretion of old 
stars, which were accreted in the same rotation direction as the 
main component. The age of the stars in the secondary compo¬ 
nent ranges from ~ 12 Gyr in the central regions down to ~ 3 
Gyr at 21", remaining nearly constant afterwards. We compare 
the mean ages of the two components within 21", which can be 
considered the upper limit for trusting age measurements in the 
main component. We And: (age)niain(^ <21") = 10 + 3 Gyr and 
(age)second(^ < 21") = 5 + 4 Gyr. The mean ages of the two 
stellar components agree within the errorbars; on the other hand, 
the radial profiles in Fig. |9] seem to indicate that the secondary 
component is systematically younger than the main galaxy by 
2 Gyr at all radial bins. Implications of this age difference, if 
real, will be discussed in Section|4] 

By comparing the mean prope rties of the stell ar popula¬ 
tions with the model predictions by iMarastonI (l2005h assuming 
a Salpeter initial mass function, we can derive the mean values 
of the stellar mass-to-light ratios of the two components in the g 
band: <M/L.)niain = 6.7 + 2.5, <M/L.)second = 3.3 + 2.5. These 
informations, combined with the luminosity of the two compo¬ 
nents measured in Sect. 13.3.21 and the adopted distance of 42.4 
Mpc, allow us to estimate their masses: Mmain = 6.2 + 2.3 ■ 10'° 
Mq and Msecond = 1.0 + 9.4 ■ IO'^Mq; the main stellar component 
is 6 times more massive than the secondary component. 

3.5. Kinematics and distribution of the Ionized gas 

The spectral decomposition routine fits also the ionized gas 
emission lines. The gas component (see Fig.|7J rotates along the 
same direction of the secondary stellar component with an am¬ 
plitude of ~ 150 km s '. The velocity dispersion decreases with 
radius: from ~ 60 km s“* in the center down to ~ 30 km s“' in 
the outskirts. The high value of velocity dispersion measured in 
the inner 1” (~ 100 km s“*) is probably due to the seeing and 
the limited spatial resolution. 

The emission lines are very weak, therefore it was necessary 
to spatially bin the data to measure their kinematics as in the 
case of the stars. For convenience, we adopted the same spatial 
binning as for the stellar component to measure the gas kine¬ 
matics. However, this “coarse” binning, prevents us to inves- 
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Fig. 7. Result of the spectral decomposition of NGC 4191. Velocity (top panels) and velocity dispersion (bottom panels) maps for the main stellar 
component (left panels), secondary counter-rotating stellar component (central panels) and counter-rotating ionized-gas component (right panels) 
are shown. Spatial bins with 5 jN lower than 25 per A are not shown for the stellar components. The horizontal bars in the color scale at the bottom 
of each panel indicate the mean error. 


NGC 4191 




Fig. 8. Line strength indices of the two counter-rotating stellar compo¬ 
nents in NGC 4191. Blue filled circles and red open diamonds refer to 
the main and secondary component, respectively. The error bars repre¬ 
sent the mean errors in the measurements. The li ne grids represent th e 
predictions of simple stellar population models dThomas et alJlTOO.^ . 
Spatial bins with S IN lower than 25 per A are not considered. 


tigate the spatial distribution of the emission lines. We there¬ 
fore re-analyzed the gas distribution on a hner spatial grid and 
measured the equivalent widths of the H/3 [O III] A4959, and 
[O III] /i5007A emission lines on this latter spatial binning. 

Figure[T0]shows the spatial distribution of the [O III] /15007A 
emission line equivalent widths, which is the most intense. Al¬ 
though the outer regions are dominated by the noise due to lower 
signal-to-noise of the data, an elongated ring-like structure (or 
maybe a spiral arm structure) with three blobs is visible within 
25". The position angle of the structure, using the two brightest 
blobs in the equivalent width map as reference, is 8.0°, which is 


very close to the position angle of the secondary component, i.e. 
7.8° (Table [T]l. This strengthens the association of the ionized- 
gas component to the secondary stellar component. 

Interestingly, a similar spiral-like structure with same size 
and orientation, is observed in the residual map of Fig.|2 consis¬ 
tently with the concentration of gas along spiral arms. Unfortu¬ 
nately, the spatial binning we adopted for the stellar population 
studies does not have enough spatial resolution to resolve the 
spiral arms and separate their stellar population from that of the 
rest of the galaxy. 

It is interesting to compare, although qualitatively, the distri¬ 
bution of emission lines of NGC 4191 with those of other stellar 
counter-rotating galaxies studied so far. The distribution of the 
ionized gas in NGC 4191 resembles those of NGC 3593 and 
NGC 5719, where the gas is also distributed along a ring-like 
structure with presence of blobs. In addition, in the case of NGC 
5719, the location of the most intense blobs coincides with that 
of the youngest stars. In the case of NGC 4550, the gas has ir¬ 
regular distribution, with no signature of ring-like or spiral-like 
structures. Also, NGC 4550 is the galaxy where the two stellar 
components are the oldest. 

4. Discussion and conciusions 

The LOSVD of NGC 4191 are complex and suggest the presence 
of two kinematically-decoupled stellar components. Their prop- 
erties have been inves tigated with the aid of a non-parametric ht 
( Fabricius et al.l 120141) and a spectral decomposition technique 
( Coccato et al.ll201 TT 

We also performed a parametric photometric decomposition 
of the surface-brightness distribution of NGC 4191. Our detailed 
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Fig. 9. Properties of the stellar populations of the counter-rotating components in NGC 4191. Left panels: median values of age (top panel), [Z/H] 
(middle panel) and [Q-/Fe] (bottom panel) within elliptical anuli for the main (red line) and secondary (blue line) stellar components. Right panels: 
two dimensional maps of the stellar population properties. The dashed lines indicate the elliptical anuli used in the computation of the median 
values presented in the left panels. Spatial bins with S /N lower than 25 per A are not considered in the analysis. 


analysis reveals the presence of three structural components, 
which were overlooked in previous studies: a Sersic bulge and 
two exponential disks of different scale lengths. Also, faint spi¬ 
ral structures (stars and dust) are visible from the residual map. 

We then linked the photometric components to the kinematic 
components unveiled by the spectroscopic analysis. The inten¬ 
sity of the double peaks in the LOSVDs allowed us to associate 
the bulge and the outer exponential disk with the main kine¬ 
matic component. Together they represent the ~ 64% of the total 
galaxy luminosity. Consequentially, we associated the inner ex¬ 
ponential disk with the secondary kinematic component, which 
represents the ~ 36%of the total galaxy light. 

The information on the photometry allowed us to constrain 
the relative light contribution of the two kinematic components 
in the spectral decomposition fit. This decreases the degeneracy 
between the various fitting parameters that characterize the two 
components (kinematics, surface brightness, and stellar popula¬ 
tions) and allows us to obtain a more reliable fit. 

We find evidences that the two kinematic components in 
NGC 4191 counter-rotate with respect to each other. The sec¬ 
ondary component rotates faster and along the same direction 
as the ionized gas. The main component rotates slower and it 
has higher velocity dispersion, consistent with hosting also bulge 
stars. 

The properties of the stellar populations of the two kinematic 
components are very similar, although the large uncertainties in 
the measurements prevent us to make strong claims. We have 
azimuthally averaged the measurements onto elliptical bins to 
reduce the statistical noise and highlight the presence of radial 
gradients. Fig.|9]gives a tentative indication that the luminosity- 
weighted ages of the stars in the secondary disk are systemati¬ 
cally lower than those of the main components at all radii. We 


can explain the age difference between the two components with 
the following scenario. 


The combination of kinematics and stellar populations gives 
us clues on the formation mechanism of NGC 4191. Our data 
favor the scenario in which a disk galaxy (represented by the 
main kinematic component, i.e., the bulge plus the outer disk) 
acquired material (either gas/stars accretion or galaxy merg¬ 
ers) from the outside onto retrograde orbit, which formed the 
counter-rotating secondary component. At first glance, this sce¬ 
nario is the sa me of those of the other counter-rotating systems 
studied so fa r (ICoccato et al.l 1201 iL 1201 3L iKatkov et al.l l201^ 
IPizzella et al.ll2014 ). However, there are some important differ¬ 
ences between NGC 4191 and these other galaxies. 


First, the age radial gradients observed in NGC 4191 are 
much steeper than those observed in other systems. More¬ 
over, the age gradient of the secondary component is negative, 
whereas in the other systems is positive. This is indicative of an 
inside-out star formation process for the secondary disk, and fa¬ 
vors the gas accretion scenario over major merger. A nearly 1:1 
merger in fact would have redistributed the stars along all direc¬ 
tions, smoothing any significant age gradient. 


Second, there are almost no differences in the metallicity and 
a-elements overabundance ratio between main and secondary 
components, whereas differences were observed in other galax¬ 
ies. In other galaxies, the scenario that better explained the dif¬ 
ferences in [Z/H] and [a/Fe] between the main and secondary 
components is gas accretion from cosmic gas reservoir or from 
a nearby gas-rich companion. The case of NGC 5719 offers 
the best example of such a scenario: we observe the on-going 
gas stripping process from it s nearby companion that forms the 
count er-rotating stellar disk (IVergani et al.l 1200^ ICoccato et al.l 
1201 ih . On the contrary, the data for NGC 4191 suggest that the 
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Fig. 10. Distribution of the [O III] 45007Aemission line. The gas is dis¬ 
tributed on a ring-like structure aligned with the photometric major axis 
of the secondary component (see text for details). 


accreted gas has origin in common with the gas that formed the 
main galaxy component. 

The observed properties of NGC 4191 can be interp reted 
in light of cosmological simulations (lAlgorrv et al.ll2014b . The 
main galaxy component is formed via accretion of gas from two 
filaments that have very similar chemical composition. As time 
passes, each filament torques the other in opposite direction, re¬ 
sulting into a “double accretion” with opposite spin. Initially, 
the accretion from the two filaments occur at the same time. The 
collisional nature of the gas ensures that only one stellar com¬ 
ponent forms, from gas whose properties is a mixture of that of 
the two filaments, and with the spin dictated by the most massive 
filament. After a while (~ 2 Gyr, considering the measured age 
difference between the two components in each radial bin), the 
accretion along one of the filaments stops. 

The accretion process continues along the counter-rotating 
filament, which wipes out the remaining gas and forms the sec¬ 
ondary counter-rotating component. The star formation occurs 
inside-out, as suggested from the negative age gradient; the new 
born stars have properties very similar to those of the main com¬ 
ponent, because the gas is a mixture of that of the second fila¬ 
ment and what remained from the first filament. Moreover, the 
formation of the secondary component must have occurred right 
after the formation of the main component and very rapidly, 
without leaving the stars of the first component enough time to 
reprocess and enrich the gas with metals and/or a elements. This 
is consistent with the relative small age difference between the 
two components and their equal and constant [Z/H] and [a/Fe] 
radial profiles. The rapid star formation process and the solar 
abundances and metallicity measured for the stellar populations, 
imply that the infalling gas was already enriched, because oth¬ 
erwise there would not be enough time to reach the measured 
values of [Z/H] and [a/Fe]. 

The hypothesis where the main component formed first is 
consistent with the observed properties of the photometric com¬ 
ponents: the position angles are slightly different and the elliptic- 
ities of the Sersic bulge and outer disk are smaller than that of the 
secondary component (inner disk, see Tab. [T]i. The most likely 
interpretation is that the bulge is triaxial and the disks have dif¬ 
ferent thickness (the external disk being the thickest). This is in 
agreement with the proposed formation scenario: the main com¬ 
ponent formed first and had more time to heat up its disk during 
the accretion of the gas that originated the secondary component. 
Gas accretion can indeed heat the host stellar di sk during the for¬ 
matio n of a counter-rotating component (i.e., iThakar & RvdenI 
119961) . 


To conclude, we would like to remark the importance of dis¬ 
entangling and measuring the properties of individual compo¬ 
nents in multi-spin galaxies to understand their formation mech¬ 
anisms. Also, the special case of NGC 4191 supports the valid¬ 
ity of the double peaked signature in the velocity dispersion field 
associated to zero-velocity rotation as selection criteria for iden¬ 
tifying counter-rotating galaxies. This is very useful for large in¬ 
tegral field spectroscopic surveys to provide a statistically com¬ 
pleted census of counter-rotating galaxies in the nearby universe. 
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Appendix A: Errors of the spectroscopic 
decomposition 

The ability of the spectral decomposition code to measure the 
kinematic and properties of the stellar populations depends on 
many parameters. First, on the differences of the two stellar com¬ 
ponents: the more separated they are in kinematics and stellar 
populations, the easier is to decouple their contribution from 
the observed galaxy spectrum. Also, the spectral resolution and 
wavelength range of the instrument and the signal-to-noise of the 
observations are important. 

We therefore simulate galaxy spectra that represent observa¬ 
tions and the instrumental set-up as close as possible. We then 
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perform the spectral decomposition on those simulated spectra 
to derive the errors on the measured quantities. 

The simulated galaxy spectra are constructed by creating and 
adding two stellar templates broadened to have 100 km s * of 
velocity dispersion. The stellar templates are both constructed 
using the same star, HD106210, observed with the VIRUS-W 
(to match the instrumental set-up), that has line strength indices 
as close as possible to the mean measured values. 

We then explore the parameters space defined by the velocity 
separation of the two component (AV), and the signal-to-noise 
ratio (S/N). For each point (AV, S/N) in the parameter space we 
construct a set of 100 simulated galaxy spectra by adding Gaus¬ 
sian noise. We then perform the spectral decomposition on those 
simulated spectra to derive the errors on the measured quantities. 
The parameter space is sampled by AV = 50,75,100,150,250 
km s ', and noise - 0.005,0.02,0.03,0.10, that correspond to 
S /N ranging from 480 to 25 per A (the mean S /N of the binned 
spectra is 90 per A) 0. 

We find that the errors on the absorption line indices does 
not vary much with velocity separation, because the stellar pop¬ 
ulations of the two components used in the simulations are the 
same (i.e. there is no degeneracy with the kinematics). There is 
also little dependence with S /N. The mean errors on the indices 
are 0.2 A 0.4 A, 0.2 A, and 0.2 A, for the HyS, Mg b, Fe5270, and 
Fe5335 indices, respectively. 

We find that the errors on velocity separation between the 
two components depend on the input velocity separation and on 
the S IN. Figure lA.ll shows the error on the recovered velocity 
separation as function of the input noise and input velocity sep¬ 
aration. For typical values measured along the galaxy kinematic 
major axis, where the velocity separation is AV > 150km s“' 
and the noise is rms < 0.05, the error on the velocity separation 
between the two components is ^ 30km s '. Also, the error on 
the recovered flux of the two components is ~ 15%. 



AV [km '] rms 



0 50 too 150 200 250 0,00 0.05 0.10 0.15 0.20 


Fig. A.l. Top figure: median value (top panel) and standard deviation 
(bottom panel) of the error on the recovered velocities as function of the 
velocity difference between the two simulated stellar components. Dif¬ 
ferent lines represent different signal-to-noise values per A: 480 (rms 
= 0.005), 120, 80, and 24 (rms = 0.1). Bottom figure: two-dimensional 
maps of the absolute velocity difference between the two stellar com¬ 
ponents (left panel) and the fit rms (right panel). 


^ We define the noise (rms) as the standard deviation of the residu¬ 
als between the galaxy observed spectrum (normalized to its median 
value) and the best fit model. The actual computation on the wave¬ 
length region used in the fit and with the robust_sigma IDL routine. 
This value of rms is converted into a proxy for the S IN per angstrom 
by 5/A = 1/rms/ VO.178, where 0.178 is the inverse dispersion in A 
pixek'. 
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